/*
 * To change this template, choose Tools | Templates
 * and open the template in the editor.
 */
package image.edge;

import image.filter.GaussianFilter;
import image.processing.ImageProcessor.CompressAlgorithm;
import image.utils.ImageUtilities;
import image.utils.ImageUtilitiesImpl;
import java.awt.image.BufferedImage;
import java.awt.image.WritableRaster;
import mask.ImageMask;

/**
 *
 * @author tommy
 */
public class HarrisCornerDetector extends EdgeDetectorImpl {
	private ImageMask horzDerivativeMask, vertDerivativeMask, gaussMask;
	private int ksize;
	
	private static double[][] horzSobelOperator = {
		  { -1d, 0d, 1d },
		  { -2d, 0d, 2d },
		  { -1d, 0d, 1d },
	};
		
	private static double[][] vertSobelOperator = {
			{-1d, -2d, -1d},
			{ 0d, 0d, 0d },
			{ 1d, 2d, 1d }
	};
		
	public HarrisCornerDetector(BufferedImage src) {
		// Applies a gaussian-filter to src before storing image
		super(new GaussianFilter(1).applyFilter(ImageUtilitiesImpl.getInstance().toGrayScale(src),
			CompressAlgorithm.LINEAR));
		this.horzDerivativeMask = new ImageMask(horzSobelOperator, 3);
		this.vertDerivativeMask = new ImageMask(vertSobelOperator, 3);
		ksize = 3;
		this.gaussMask = new ImageMask(buildGaussianKernel(),ksize);
		
		// resets image "ret" and copies src into it, with 24bit-rgb format
		ret = ImageUtilitiesImpl.getInstance().toGrayScale(src);
		ret = ImageUtilitiesImpl.getInstance().toRGB(ret);
	}
	
	public BufferedImage FindEdges() {
		Double[] iX, iY, a, b, c; // image directional derivatives
		double[] iX2, iY2, iXiY; // image squared and cross directional derivatives
		double[] hXY; // cornerness H(x,y)
		int w, h, wh;
		
		BufferedImage A, B, C; // Ix2, Iy2 and IxIy images
		WritableRaster wA, wB, wC, outRaster;
		
		w = ret.getWidth();
		h = ret.getHeight();
		wh = w*h;
		
		// Applies Sobel Operator to compute finite differentiation in 'X' and 'Y'
		// does not process padded border edge
		iX = horzDerivativeMask.applyMask(src);
		iY = vertDerivativeMask.applyMask(src);
		
		iX2 = new double[iX.length];
		iY2 = new double[iY.length];
		iXiY = new double[iX.length];
		hXY = new double[iX.length];
		
		// squared and cross derivatives
		for (int i = 0 ; i < iX2.length ; i++) {
			iX2[i] = iX[i] * iX[i];
			iY2[i] = iY[i] * iY[i];
			iXiY[i] = iX[i] * iY[i];
		}
		ImageUtilities imgUtil = ImageUtilitiesImpl.getInstance();
		
		A = imgUtil.addPaddingEdge(new BufferedImage(w, h, BufferedImage.TYPE_BYTE_GRAY), 1);
		B = imgUtil.addPaddingEdge(new BufferedImage(w, h, BufferedImage.TYPE_BYTE_GRAY), 1);
		C = imgUtil.addPaddingEdge(new BufferedImage(w, h, BufferedImage.TYPE_BYTE_GRAY), 1);
		
		wA = A.getRaster();
		wB = B.getRaster();
		wC = C.getRaster();
		outRaster = ret.getRaster();
		
		wA.setSamples(1, 1, w, h, 0, iX2);
		wB.setSamples(1, 1, w, h, 0, iXiY);
		wC.setSamples(1, 1, w, h, 0, iY2);
		
		/*gaussFilter = new GaussianFilter(1);
		A = gaussFilter.applyFilter(A, CompressAlgorithm.LINEAR);
		B = gaussFilter.applyFilter(B, CompressAlgorithm.LINEAR);
		C = gaussFilter.applyFilter(C, CompressAlgorithm.LINEAR);
		*/
		
		a = gaussMask.applyMask(A);
		b = gaussMask.applyMask(B);
		c = gaussMask.applyMask(C);
		
		
		// compute cornerness H(x,y)
		for (int i = 0 ; i < wh ; i++) {
			//hXY[i] = a[i]*c[i] - b[i]*b[i] - 0.04*(Math.pow(a[i]+c[i], 2));
			hXY[i] = a[i]*c[i] - b[i]*b[i] / (a[i] + c[i] + Double.MAX_VALUE);
		}
		
		findLocalMaxima(hXY, w, h);
		
		for (int i = 0 ; i < hXY.length ; i++) {
			if (hXY[i] > 0) {
				outRaster.setSample(i%w, i/w, 0, hXY[i]);
			}
		}
		
		
		return ret;
	}
	
	private void findLocalMaxima(double[] hXY, int w, int h) {
		
		for (int i = 0 ; i < w * h ; i++) {
			// check if index targets at image borders
			if (i < w || i >= w*h - w || 
				i % w == 0 || i % w == w - 1) {
				hXY[i] = 0d;
			} else {
				
				if (!isCornerLocalMax(hXY, i, w, h)) {
					hXY[i] = 0;
				} else {
					hXY[i] = 255;
				}
			}
		}
	}
	
	// idx: flat index
	private boolean isCornerLocalMax(double[] hXY, int idx, int w, int h) {
		boolean isMax = true;
		int voffset; // vertical offset to move to target x,y location in 3x3 kernel
		int hoffset; // horizontal offset
		
		for (int j = 0 ; j < 3 && isMax == true; j++) {
			voffset = (j - 1) * w;
			for (int i = 0 ; i < 3 && isMax == true ; i++) {
				hoffset = i - 1;
				
				if (i != 1 || j != 1) {
					if (hXY[idx] <= hXY[idx + voffset + hoffset]) {
						isMax = false;
					}
				}
			}
		}
		
		return isMax;
	}
	
	private double[][] buildGaussianKernel() {
		double[][] gaussianKernel;
		double xmapIndex, ymapIndex; // transform (x,y) to gaussian kernel coords
		double sigma2, exp, gaussValue, sum;
		
		sum = 0;
			
		gaussianKernel = new double[ksize][ksize];
		sigma2 = 1; //^2
		
		for (int y = 0 ; y < ksize ; y++) {
			for (int x = 0 ; x < ksize ; x++) {
				xmapIndex = x - ksize/2;
				ymapIndex = y - ksize/2;
				exp = Math.exp(-((xmapIndex*xmapIndex + ymapIndex*ymapIndex) / (2 * sigma2) ));
				gaussValue = (1d / ((2 * Math.PI) * sigma2)) * exp;
				gaussianKernel[y][x] = gaussValue;
				sum += gaussValue;
			}
		}
		
		// Operation to make kernel elements sum to 0
		sum /= (ksize*ksize);
		
		for (int y = 0 ; y < ksize ; y++) {
			for (int x = 0 ; x < ksize ; x++) {
				gaussianKernel[y][x] -= sum;
			}
		}
		return gaussianKernel;
	}	
}
